#include "../../headers/structsystem.h"

using namespace std;

struct Point
{
    int id;
    double x, y, z;
};

struct Element
{
    int id;
    int mat;
    int typ;
    int sec;
    int p1, p2;
};

void importNodes(const string &fname, vector<Point> &v1)
{
    ifstream fin;
    fin.open(fname.c_str());
    if (!fin.is_open())
    {
        cerr << fname << ": open file failed" << endl;
        exit(-1);
    }

    string line;
    getline(fin, line);

    char comma;
    while (true)
    {
        Point p;
        fin >> p.id >> comma >> p.x >> comma >> p.y >> comma >> p.z;
        if(!fin.good()) break;
        v1.push_back(p);
    }
    fin.close();
}

void importElements(const string &fname, vector<Element> &v1)
{
    ifstream fin;
    fin.open(fname.c_str());
    if (!fin.is_open())
    {
        cerr << fname << ": open file failed" << endl;
        exit(-1);
    }

    string line;
    getline(fin, line);

    char comma;
    while (true)
    {
        Element e;
        fin >> e.id >> comma >> e.mat >> comma >> e.typ >> comma >> e.sec
            >> comma >> e.p1 >> comma >> e.p2;
        if(!fin.good()) break;
        v1.push_back(e);
    }
    fin.close();
}

void saveExternalForce(const std::string &path, Particle* p)
{
    // create file
    string fname = path + "/particle_" + to_string(p->id())+ ".csv";
    fstream fout;
    fout.open(fname, ios::out | ios::app);

    // write header
    string header = "PID, Fx, Fy, Fz, Mx, My, Mz";
    fout << header << endl;

    // write data
    p->outputReactionForce(fout);

    // close file
    fout.close();
}

int main()
{
    // create system
    StructSystem system = StructSystem(1);
    system.setJobname("Star");

    // solve parameters
    double endTime = 100.0;
    double zeta = 1;
    double d = 6.0;
    double rtime = 50.0;

    // import structure parameters
    string fnode = system.workdir() + "/" + "nodes.csv";
    vector<Point> points;
    importNodes(fnode, points);

    string felems = system.workdir() + "/" + "elements.csv";
    vector<Element> elements;
    importElements(felems, elements);

    // create particles
    map<int, Particle*> particles;
    for (int i = 0; i < points.size(); i++)
    {
        double x = points[i].x;
        double y = points[i].y;
        double z = points[i].z;
        Particle* p = new Particle(points[i].id, x, y, z);
        system.addParticle(p);
        particles[points[i].id] = p;
    }

    // create material
    LinearElastic mat = LinearElastic(1);
    mat.setE(1.0E6);
    mat.setDensity(20);

    // create section
    CustomSectionParas paras = {.A=0.1, .Sy=0, .Sz=0, .SHy=0, .SHz=0,
                          .CGy=0, .CGz=0.0, .Iyy=0.01, .Izz=0.01, .Iyz=10};
    CustomSection sect = CustomSection(1, paras);

    // create elements
    for (int i = 0; i < elements.size(); i++)
    {
        int id = elements[i].id;
        int p1 = elements[i].p1;
        int p2 = elements[i].p2;
        Link3DLD* e = new Link3DLD(id, particles[p1], particles[p2],
                           &mat, &sect);
        system.addElement(e);
    }

    // constraints
    DOF v1 = {.key={true, true, true, true, true, true},
              .val={0, 0, 0, 0, 0, 0}};
    system.addConstraint(8, v1);
    system.addConstraint(9, v1);
    system.addConstraint(10, v1);
    system.addConstraint(11, v1);
    system.addConstraint(12, v1);
    system.addConstraint(13, v1);
    system.info();

    string path = system.workdir() + "/" + system.jobname() + "/model";
    system.saveModel(path);
    StdArray6d d1 {};

    // setting parameters
    double h = system.autoTimeStep();
    system.setDampCoeff(zeta);
    int nStep = ceil(endTime / h);

    int interval = ceil(0.1 / h);
    cout << "##### start calculating ######" << endl;
    for (int i = 0; i <= nStep; i++)
    {
        d1[2] = -d * (double)i / nStep;
        if (i == 0)
        {
            // system.addExternalForce(1, d1);
            system.setParticleDisplace(1, d1);
            // particles[1]->setDisplace(d1);
            system.solve(h, zeta, true);
            system.setInternalForce();
        }
        else
        {
            system.solve(h, zeta, false);
            system.clearParticleForce();
            // system.addExternalForce(1, d1);
            // particles[1]->setDisplace(d1);
            system.setParticleDisplace(1, d1);
            system.setInternalForce();
        }
        // save results
        if (i % interval == 0)
        {
            string path = system.workdir() + "/" + system.jobname() + "/" +
                          to_string(i * h);
            system.saveParticleResult(path);
            system.saveElementResult(path);
            system.saveSupportReact(path);
            saveExternalForce(path, particles[1]);
            saveExternalForce(path, particles[2]);
        }
    }
    cout << "##### finished ######" << endl;

    system.releaseContainers();
    return 0;
}